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We present a comprehensive study of Vicsek-style self-propelled particle models in two and three 
space dimensions. The onset of collective motion in such stochastic models with only local alignment 
interactions is studied in detail and shown to be discontinuous (first-order like). The properties of the 
ordered, collectively moving phase are investigated. In a large domain of parameter space including 
the transition region, well-defined high-density and high-order propagating solitary structures are 
shown to dominate the dynamics. Far enough from the transition region, on the other hand, these 
objects are not present. A statistically-homogeneous ordered phase is then observed, which is 
characterized by anomalously-strong density fluctuations, superdiffusion, and strong intermittency. 

PACS numbers: 64.60Cn, 05.70Ln 



I. INTRODUCTION 

Collective motion phenomena in nature have attracted 
the interest of scientists and other authors for quite a 
long time The question of the advantage of hving 
and moving in groups, for instance, is a favorite one 
among evolutionary biologists . In a different perspec- 
tive, physicists are mostly concerned with the mecha- 
nisms at the origin of collective motion, especially when 
it manifests itself as a true, non-trivial, emerging phe- 
nomenon, i.e. in the absence of some obvious cause like 
the existence of a leader followed by the group, a strong 
geometrical constraint forcing the displacement, or some 
external field or gradient felt by the whole population. 
Moreover, the ubiquity of the phenomenon at all scales, 
from intra-cellular molecular cooperative motion to the 
displacement in group of large animals, raises, for physi- 
cists at least, the question of the existence of some uni- 
versal features possibly shared among many different sit- 
uations. 

One way of approaching these problems is to construct 
and study minimal models of collective motion: if uni- 
versal properties of collective motion do exist, then they 
should appear clearly within such models and thus could 
be efficiently determined there, before being tested for 
in more elaborate models and real-world experiments or 
observations. Such is the underlying motivation of re- 
cent studies of collective motion by a string of physi- 
cists d, 0, H, H, 0, Hi ■ Among them, the group of Tamas 
Vicsek has put forward what is probably the simplest pos- 
sible model exhibiting collective motion in a non-trivial 
manner. 

In the "Vicsek model" [9], point particles move off- 



lattice at constant speed vq, adjusting their direction of 
motion to that of the average velocity of their neighbors, 
up to some noise term accounting for external or internal 
perturbations (see below for a precise definition). For a 
finite density of particles in a finite box, perfect align- 
ment is reached easily in the absence of noise: in this 
fluctuation-less collective motion, the macroscopic veloc- 
ity equals the microscopic one. On the other hand, for 
strong noise particles are essentially non-interacting ran- 
dom walkers and their macroscopic velocity is zero, up 
to statistical fluctuations. 

Vicsek et al. showed that the onset of collective motion 
occurs at a finite noise level. In other words, there exists, 
in the asymptotic limit, a fluctuating phase where the 
macroscopic velocity of the total population is, on aver- 
age, finite. Working mostly in two space dimensions, they 
concluded, on the basis of numerical 0, Q simulations, 
that the onset of this ordered motion is well described 
as a novel non-equilibrium continuous phase transition 
leading to long range order, at odds with equilibrium 
where the continuous XY symmetry cannot be sponta- 
neously broken in two space dimensions and below p^ . 
This brought support to the idea of universal properties 
of collective motion since the scaling exponents and func- 
tions associated to such phase transitions are expected to 
bear some degree of universality, even out of equilibrium. 

The above results caused a well-deserved stir and 
promp_ted_ a_ large number of studies at various lev- 
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In particular, two of us showed that the onset of collective 
motion is in fact discontinuous [s^, and that the origi- 
nal conclusion of Vicsek et al. was based on numerical 
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results obtained at too small sizes [1, Q . More recently, 
the discontinuous character of the transition was chal- 
lenged in two publications, one by the Vicsek group [3§| 
and another by Aldana et al. [40| . 

Here, after a definition of the models involved (Sec- 
tion|TT|, we come back, in Section|TTTl to this central issue 
and present a rather comprehensive study of the onset of 
collective motion in Vicsek-style models. In Section IIVI 
we describe the ordered, collective motion phase. Sec- 
tion |V] is devoted to a general discussion of our results 
together with some perspectives. Most of the numeri- 
cal results shown were obtained in two space dimensions, 
but we also present three-dimensional results. Wherever 
no explicit mention is made, the default space dimension 
is two. Similarly, the default boundary conditions are 
periodic in a square or cubic domain. 



when particles try to follow the locally-averaged direc- 
tion of motion. One could argue on the other hand that 
most of the randomness stems from the evaluation of each 
interaction between particle i and one of its neighbors, 
because, e.g., of perception errors or turbulent fluctua- 
tions in the medium. This suggests to replace Eq.([T]) by: 



(3) 



where ^ is a random unit vector and A/i is the number of 
particles in Si. It is easy to realize that this "vectorial" 
noise acts differently on the system. While the intensity 
of angular noise is independent from the degree of local 
alignment, the influence of the vectorial noise decreases 
with increasing local order. 



II. THE MODELS 

A. Vicsek model: angular noise 

Let us first recall the dynamical rule defining the Vic- 
sek model (eJ] . Point particles labeled by an integer index 
i move off-lattice in a space of dimension d with a veloc- 
ity Vi of fixed modulus vq = \vi\. The direction of motion 
of particle i depends on the average velocity of all par- 
ticles (including i) in the spherical neighbourhood Si of 
radius centered on i. The discrete-time dynamics is 
synchronous: the direction of motion and the position of 
all particles are updated at each timestep At, in a driven, 
overdamped manner: 



vlit + At) =wo (7^,,o^9) 



where i? is a normalization operator 



(1) 



J|) and 



C. Repulsive force 

In the original formulation of the Vicsek model as well 
as in the two variants defined above, the only interaction 
is alignment. In a separate work [22| . we introduced a 
two-body repulsion/attraction force, to account for the 
possibility of maintaining the cohesion of a flock in an 
infinite space (something the Vicsek model does not al- 
low). Here, we only study models without cohesion. Nev- 
ertheless, we have considered, in the following, the case 
of a pairwise repulsion "force" , to estimate in particular 
the possible infiuence of the absence of volume exclusion 
effects in the basic model, which leaves the local den- 
sity actually unbounded. We thus introduce the short 
ranged, purely repulsive interaction exerted by particle j 
on particle i: 



-Cij x[l + exp(|fj - fi\/rc - 2)]" 



(4) 



TZr^ performs a random rotation uniformely distributed where e*y is the unit vector pointing from particle i to j 



around the argument vector: in d = 2, TZ^fV is uniformely 
distributed around H inside an arc of amplitude 27r 77; in 
d = 3, it lies in the solid angle subtended by a spherical 
cap of amplitude 47r 77 and centered around w. The par- 
ticles positions fl are then simply updated by streaming 
along the chosen direction as in 



and Tc < tq is the typical repulsion range. Equations ([T]) 
and ([3]) are then respectively generalized to 



Vi{t + At) = vq (7^^o^9) 



(5) 



fl{t + At) = nit) + At vl{t + At) 



(2) 



and 



Note that the original updating scheme proposed by Vic- 
sek et al. in Q defined the speed as a backward differ- 
ence, although we are using a forward difference. The 
simpler updating above, now adopted in most studies of 
Vicsek-style models, is not expected to yield different re- 
sults in the asymptotic limit of infinite size and time. 



Vi{t + At) ^VQ-d 



Y.v-it) + pJ2k+'^^^^ 



(6) 



where (3 measures the relative strength of repulsion with 
respect to alignment and noise strength. 



A different noise term: vectorial noise 



D. Control and order parameters 



The "angular" noise term in the model defined above 
can be thought of as arising from the errors committed 



The natural order parameter for our polar particles is 
simply the macroscopic mean velocity, conveniently nor- 
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FIG. 1: (color online) Typical behavior across the on- 
set of collective motion for moderate size models (p = 2, 
no = 0.5, L = 64) with angular noise (black circles) and 
vectorial noise (red triangles). The reduced noise amplitude 
£ = l — r]/rjt is shown in abscissas (transition points estimated 
at Tjt = 0.6144(2) - vectorial noise - and r]t = 0.478(5) - an- 
gular noise), (a) Time-averaged order parameter {(p{t))t. (b) 
Binder cumulant G. (c) Variance a of (p; (d) Order parame- 
ter distribution function P at the transition point. Bimodal 
distribution for vectorial noise dynamics (red dashed line), 
unimodal shape for angular noise (black solid line). Time 
averages have been computed over 3 • 10^ timesteps. 



malized by the microscopic velocity vq 



(7) 



where stands for the average over the whole popula- 
tion. Here, we mostly consider its modulus ip{t) = 11^(^)1, 
the scalar order parameter. 

In the following, we set, without loss of generality. 
At — 1 and rg = 1, and express all time- and length- 
scales in terms of these units. Moreover repulsive force 
will be studied by fixing Tc — 0.127 and f3 = 2.5. 

This leaves us with two main parameters for these 
models: the noise amplitude rj and the global density 
of particles p. Recently, the microscopic velocity vq has 
been argued to play a major role as well [S^. All three 
parameters (?], p, and vq) are considered below. 



III. ORDER-DISORDER TRANSITION AT THE 
ONSET OF COLLECTIVE MOTION 

As mentioned above, the original Vicsek model at- 
tracted a lot of attention mostly because of the conclu- 




FIG. 2: (color online) FSS analysis of angular noise dynam- 
ics (p — 2, vo = 0.5, time averages computed over 2 • 10^ 
timesteps). Time-averaged order parameter (a) and Binder 
cumulant (b) as a function of noise for various system sizes 
L. (c) Piece of an order parameter time series close to the 
transition point (L = 256, 77 = 0.476). 



sions drawn from the early numerical studies [5I, ^ : the 
onset of collective motion was found to be a novel contin- 
uous phase transition spontaneously breakiiig rotational 
symmetry. However, it was later shown in t38| that be- 
yond the typical sizes considered originally the discon- 
tinuous nature of the transition emerges, irrespective of 
the form of the noise term. Recently, the discontinuous 
character of the transition was argued to disappear in the 
limit of small vq ^SOj . We now address the problem of the 
nature of the transition in full detail. 

Even though there is no rigorous theory for finite- 
size-scaling (FSS) for out-of-equilibrium phase transi- 
tions, there exists now ample evidence that one can 
safely rely on the knowledge gained in equilibrium sys- 
tems m m, H. The FSS approach ^ii, '45^ involves 
the numerical estimation of various moments of the or- 
der parameter distribution as the linear system size L 
is systematically varied. Of particular interest are the 
variance 

and the so-called Binder cumulant 



2\2 ' 



(8) 



where (•)( indicates time average. The Binder cumulant 
is especially useful in the case of continuous phase tran- 
sitions, because it is one of the simplest ratio of moments 
which takes a universal value at the critical point r^t, 
where all the curves G{r],L), obtained at different sys- 
tem sizes L, cross each other. At a first-order transition 
point, on the other hand, the Binder cumulant exhibits 
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FIG. 3: (color online) Transition to collective motion with 
short range repulsive interactions. Left panels: angular noise. 
Right panels: vectorial noise. (a,b) : order parameter vs noise 
amplitude at different system sizes. (c,d): Binder cumulant 
G as a function of noise amplitude, (p = 2, vq = 0.3, time 
averages carried over 10^ timesteps). 
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FIG. 4: (color online) Transition to collective motion in three 
spatial dimensions. Left panels: angular noise. Right pan- 
els: vectorial noise. (a,b): time-averaged order parameter vs. 
noise amplitude at different system sizes. (c,d): Binder cu- 
mulant G as a function of noise amplitude, (p = 0.5, vq = 0.5, 
time averages carried over 10^ timesteps). 



a sharp drop towards negative values This mini- 

mum is due to the simuhaneous contributions of the two 
phases coexisting at threshold. Moreover, it is easy to 
compute that G{ri,L) « 2/3 in the ordered phase, while 
for a disordered state with a continuous rotational sim- 
metry one has G{ri, L) « 1/3 in d = 2 and G(?7, L) w 4/9 
in d = 3. 



A. Overture 

As an overture, we analyze systems of moderate size 
in two dimensions {N w 10^ particles) at the density 
/o = 2, typical of the initial studies by Vicsek et ai, but 
with the slightly modified update rule ^ and for both 
angular and vectorial noise. The microscopic velocity is 
set to vq = 0.5. 

For angular noise, the transition looks indeed contin- 
uous, as found by Vicsek et at. On the other hand, 
the time averaged scalar order parameter {(p)t displays 
a sharp drop for vectorial noise, and the Binder cumu- 
lant exibits a minimum at the transition point, indicating 
a discontinuous phase transition (Fig. [T^-b). Simultane- 
ously, the variance is almost delta-peaked. The difference 
betwen the two cases is also recorded in the probability 
distribution function (PDF) of ip which is bimodal (phase 
coexistence) in the vectorial noise case (Fig. [T}:-d) . 

The qualitative difference observed upon changing the 



way noise is implemented in the dynamics is, however, 
only a finite-size effect. As shown in the transition 
in the angular noise case reveals its asymptotic discon- 
tinuous character provided large-enough system sizes L 
are considered (Fig.[2^-b). Remaining for now at a qual- 
itative level, we show in Fig. [2t a typical time-series of 
the order parameter for the angular noise case in a large 
system in the transition region. The sudden jumps from 
the disordered phase to the ordered one and vice-versa 
are evidence for metastability and phase coexistence. 

Note that the system size beyond which the transi- 
tion reveals its discontinuous character for the angular 
noise case at density p = 2 and velocity vq = 0.5 — the 
conditions of the original papers by Vicsek et al. — is of 
the order of L = 128, the maximum size then consid- 
ered. It is clear also from Fig. [T] that the discontinuous 
nature of the transition appears earlier, when increasing 
system size, for vectorial noise than for angular noise. 
Thus, finite-size effects are stronger for angular noise. 
The same is true when one is in presence of the repulsive 
interactions (Fig[3]). Finally, the same scenario holds in 
three space dimensions, with a discontinous phase tran- 
sition separating the ordered from the disordered phases 
for both angular and vectorial noise (Fig. [4]). 

Before proceeding to a study of the complete phase 
diagram, we detail now how a comprehensive FSS study 
can be performed on a particular case. 
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FIG. 5: Correlation time r of the order parameter near the 
transition point for vectorial noise dynamics with repulsion. 
System parameters are p = 2, vq = 0.5 and rj rjt. (a) Time 
correlation function C{t) at L = 128. The lin-log inset shows 
the exponential decay, (b) Correlation time r as a function 
of system size L. The dashed line marks linear growth with 
L. Correlation functions were computed on samples of ~ lO'' 
realizations for typically 10'' timesteps. 



B. Complete FSS analysis 

For historical reasons, the following study has been 
performed on the model with vectorial noise and repul- 
sive force (Eq. [S]) . It has not been repeated in the simpler 
case of the "pure" Vicsek model because its already high 
numerical cost would have been prohibitive due to the 
strong finite size effects. 

As a first step, we estimated the correlation time t{L), 
whose knowledge is needed to control the quality of time- 
averaging: the duration T of numerical simulations has 
been taken much larger than t{L) (T = lOOr in the 
largest systems, but typically lOOOOr for smaller sizes). 
Moreover, r is also useful to correctly estimate the sta- 
tistical errors on the various moments (as {ip)t, cr, and G 
) of the PDF of the order parameter, for which we used 
the Jackknife procedure [47[. The correlation time was 
estimated near the transition (where it is expected to be 
largest) as function of system size L measuring the expo- 
nential decay rate of the correlation function (Fig. [5^) 



C{t) ^ {cpitoMto +t))t - {ip{to))l 



exp 



t 



(9) 



We found t to vary roughly linearly with L (see Fig.tSjj). 
It is interesting to observe that at equilibrium, one would 
expect r to scale as 

T = Li expinL'^-^) 

where k is the surface tension of the metastable state. 
Therefore, our result implies a very small or vanishing 
surface tension k <C a situation reminiscent of 

observations made in the cohesive case [22], where the 
surface tension of a cohesive droplet was found to vanish 
near the onset. 

Following Borgs and Kotecky [i^, the asymptotic co- 
existence point rjt (i.e. the first order transition point) 




128 256 



FIG. 6: (color online) FSS analysis of vectorial dynamics with 
short range repulsive force (p = 2, vo — 0.3). Convergence of 
the finite-size transition points measured from different mo- 
ments of the order parameter FSS to the asymptotic transi- 
tion point rjt (see fig. |3]). Upper panels: finite size transition 
points estimated from (a) time average, (b) variance and (c) 
Binder cumulant. The horizontal dashed line marks the es- 
timated asymptotic threshold rjt = 0.5544(1). Lower panels: 
scaling of the finite size reduced noise e = 1 — rj/rjt transition 
point, (d) Exponential convergence for the jump location in 
the time-averaged order parameter, (e) Power-law behavior 
of the variance peak position, (f) Power-law behavior of the 
Binder cumulant minimum. The dashed lines in (e-f) mark 
the estimated exponents 7<t = 7g = 2. 



can be determined from the asymptotic convergence of 
various moments of the order parameter PDF. First, the 
observed discontinuity in {(p(t))t, located at rj^[L), is ex- 
pected to converge exponentially to rjt with L. Second, 
the location of the susceptibility peak ria-(L) — which 
is the same as the peak in a provided some fluctuation- 
dissipation relation holds, see the Appendix — also con- 
verges to rjt, albeit algebraically with an exponent j^r- 
Third, the location of the minimum of G, TjciL), is also 
expected to converge algebraically to rjt with with an ex- 
ponent 7g = 7^. 

Interestingly, the value taken by these exponents ac- 
tually depends on the number of phases and of the di- 
mension d of the system: for two-phase coexistence one 
has 7g = 7(7 = 2c?, while for more than two phases 
7g = 7o- = d. In Figure O we show that our data are 
in good agreement with all these predictions. The three 
estimates of rjt are consistent with each other within nu- 
merical accuracy. Moreover ri^{L) is found to converge 
exponentially to the transitional noise amplitude, while 
both T]^(L) and t]g{L) show algebraic convergence with 
an exponent close to 2. This agrees with the fact that 
due to the continuous rotational symmetry, the ordered 
phase is degenerate and amounts to an infinite number 
of possible phases. 
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FIG. 7: (color online) Hysteresis in three spatial dimensions 
with vectorial noise, (a) order parameter vs noise strength 
along the hysteresis loop observed with a ramp rate of 2 • 
10"'' per time step (p = 1/2, vo = 0.5, L = 32). Empty 
circles mark the path along the adiabatic increase of noise 
amplitude, full triangles for adiabatic decrease, (b) nucleation 
times from the disordered phase to the ordered phase (r| , left 
curves) and vice- versa (r|, right curves) for two system sizes 
(other parameters as in (a)). Each point is averaged over 1000 
realizations. 



C. Hysteresis 

One of the classical hallmarks of discontinuous phase 
transitions is the presence, near the transition, of the 
hysteresis phenomenon: ramping the control parameter 
at a fixed (slow) rate up and down through the transi- 
tion point, a hysteresis loop is formed, inside which phase 
coexistence is manifest (see Fig. [7^ for the d = 3 case 
with vectorial noise). The size of such hysteresis loops 
varies with the ramping rate. An intrinsic way of as- 
sessing phase coexistence and hysteresis is to study sys- 
tematically the nucleation time r| needed to jump from 
the disordered phase to the ordered one, as well as tj^ 
the decay time after which the ordered phase falls into 
the disordered one. Fig. [7|d shows, in three space dimen- 
sions, how these nucleation and decay vary with rj at two 
different sizes. A sharp divergence is observed, corre- 
sponding to the transition point. At a given time value 
T, one can read, from the distance between the "up" and 
the "down" curve, the average size of hysteresis loops for 
ramping rates of the order of 1/r. 



D. Phase diagram 

The above detailed FSS study would be very tedious to 
realize when varying systematically the main parameters 
rj, p, and vq, as well as the nature of the noise and the 
presence or not of repulsive interactions. From now on, 
to characterize the discontinuous nature of the transition, 
we rely mainly on the presence, at large-enough system 
sizes L, of a minimum in the variation of the Binder 
cumulant G with rj (all other parameters being fixed). 
We call L* the crossover size marking the emergence of 
a minimum of G(r]). 
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FIG. 8: (color online) Asymptotic phase diagrams for the 
transition to collective motion, (a) Two space dimensions: 
threshold amplitude r\t for angular noise as a function of 
density p at iiq = 0.5. Inset: Log-log plot to compare the 
low density behavior with the mean field predicted behavior 
r\t ~ yp (dashed red line), (b) As in panel (a), buth with 
vectorial noise dynamics, (c) Noise-density phase diagram in 
three dimensions for vectorial noise dynamics at fixed veloc- 
ity Vq — 0.5. In the log-log inset the transition line can be 
compared with the predicted behavior rjt ~ p^^^ (dashed red 
line), (d) Two space dimensions: threshold amplitude % for 
angular noise as a function of particle velocity vq at fixed den- 
sity p = 1/2 (black circles) and p = 1/8 (red triangles). The 
horizontal dashed line marks the noise amplitude considered 
in Ref. [39| (see Section HIlFl l . 



We are now in the position to sketch the phase diagram 
in the {rj, p, vq) parameter space. The numerical protocol 
used is, at given parameter values, to run a large-enough 
system so that the discontinuous character of the transi- 
tion is seen (i.e. L > L*). For larger sizes, the location of 
the transition point typically varies very little, so that for 
most practical purposes, locating the (asymptotic) tran- 
sition point from systems sizes around L* is satisfactory. 

The results presented below are in agreement with sim- 
ple mean- field-like arguments in the diluted limit: in the 
small-p regime, one typically expects that the lower the 
density, the lower the transitional noise amplitude rjt- 
Indeed, for Atwo of the order of or not much smaller 
than the interaction range rg and in the low-density limit 
p <C I/t'o'', the system can be seen as a dilute gas in which 
particles interact by short range ordering forces only. In 
this regime, the persistence length of an isolated parti- 
cle (i.e. the distance travelled before its velocity loses 
correlation with its initial direction of motion) varies like 
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FIG. 9: (color online) First-order transition for angular noise 
dynamics at high (left panels) and low velocity no (right pan- 
els). Typical averaging time is « 10^ timesteps. (a) Time- 
averaged order parameter and (c) Binder cumulant at large 
particle velocity for angular noise in two spatial dimensions at 
increasing velocities and L > L* (p = 2). (b) Time averaged 
order parameter and (d) Binder cumulant for vo = 0.05 and 
two increasing system sizes (p = 1/2). 



vo/rj. To allow for an ordered state, the noise amplitude 
should be small enough so that the persistence length re- 
mains larger than the average inter-particle distance, i.e. 
Xj p^l'^. Thus the transition noise amplitude is expected 
to behave as 



»7t 



(10) 



In 0, it was indeed found that r\i ^ p" with a ~ 
\ in two dimensions. Our own data (Fig. [SK-c) now 
confirm Eq. (fTO|) for both the angular and vectorial noise 
in two and three spatial dimensions, down to very small 
p values. The data deviate from the square-root behavior 
as the average inter-particle distance gets of the order of 
or smaller than the interaction range. 

Finally, we also investigated the transition line when uq 
is varied (Fig. [8ji). For the vectorial noise case, at fixed 
density, the threshold noise value ryt is almost constant 
(data obtained at p = i, not shown). For the angular 
noise, in the small wq limit where the above mean-field 
argument does not apply, we confirm the first-order char- 
acter of the phase transition down to ~ 0.05 for both 
angular and vectorial noise (Fig. [9b, d). For even smaller 
values of wq, the investigation becomes numerically too 
costly (see section below). Note that 771 seems to be finite 
when ^ 0"*", a limit corresponding to the XY model on 
a randomly connected graph. Still for angular noise, the 
large velocity limit is also difficult to study numerically. 
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FIG. 10: (color online) Crossover system size L* above which 
the discontinuous character of the transition appears (as testi- 
fied by the existence of a minimum in the G{r\) curve). Black 
circles: angular noise. Red triangles: vectorial noise, (a) L* 
vs p for = 0.5; Inset: the low density behavior in log-log 
scales; the dashed line marks a power-law divergence propor- 
tional to 1/y^. (b) L* vs vq at fixed density (open symbols: 
p = 1/2; filled symbols p = 2.0). 



Again, we observe that the transition is discontinuous as 
far as we can probe it, i.e. = 20 (Fig. [n^,b). 



E. Special limits and strength of finite-size effects 

We now discuss particular limits of the models above 
together with the relative importance of finite-size effects. 
Recall that these are quantified by the estimated value 
of the crossover size L* beyond which the transition ap- 
pears discontinuous. All the following results have been 
obtained for d = 2. Partial results in three dimensions in- 
dicate that the same conclusions should hold there. Keep 
in mind that in all cases reported, the transition is dis- 
continuous. We are just interested here in how large a 
system one should use in order to reach the asymptotic 
regime. 

Fig. [TUb shows that finite-size effects are stronger for 
angular noise than for vectorial noise for all densities p at 
which we are able to perform these measurements. Note 
in particular that at p = 2, the density originally used 
by Vicsek el ai, L* ~ 128 for angular noise, while it is 
very small for vectorial noise, confirming the observation 
made in Section fill AI 

In the small-p limit, the discontinuous character of the 
transition appears later and later, with L* roughly di- 
verging as 1/ (inset of Fig. [TOk). Note that this means 
that in the small-p limit one needs approximately the 
same number of particles to start observing the disconti- 
nuity. 

The large-p limit reveals a difference between angular 
and vectorial noise: while L* remains small for vectorial 
noise, it seems to diverge for angular noise (Fig. [TUk). 
making this case difficult to study numerically. 

We also explored the role the microscopic velocity vq 
in the strength of finite-size effects. Qualitatively, the 
effects observed are similar to those just reported when 
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the density is varied (Fig. [TOb). In the small vq limit, 
we record a strong increase oi L* as vq ^ for both 
types of noise. In the large velocity limit, L* decreases 
for vectorial noise, whereas it increases for angular noise. 



F. Summary and discussion 

The summary of the above lengthy study of the order- 
disorder transition in Vicsck-like models is simple: for 
any finite density p, any finite velocity Wq, and for both 
types of noise introduced, the transition is discontinuous. 
This was observed even in the numerically-difficult limits 
of large or small p or vq- These results contradict recent 
claims made about the angular noise case (original Vicsek 
model). We now comment on these claims. 

Vicsek and co-workers [s^, showed that keeping the 
density and the noise intensity fixed, a qualitative change 
is observed when decreasing vq: for not too small vq val- 
ues, in the ordered phase, particles diffuse anisotropically 
(and the transition is discontinuous), while diffusion be- 
comes isotropic at small vq, something interpreted as a 
sign of a continuous transition in this region. Rather 
than the convoluted arguments presented there, what 
happens is in fact rather simple: decreasing vq at fixed 
p and rj, one can in fact cross the transition line, pass- 
ing from the ordered phase (where particles obviously 
diffuse anisotropically due to the transverse superdiffu- 
sive effects discussed in Section IIV Cp to the disordered 
phase. Our Figure [Sji, obtained in the same conditions 
as in [39j (apart from harmless change of time updating 
rule), shows that if one keeps 77 = 0.1 (as in [3^), one 
crosses the transition line at about vq ~ 0.1, the value 
invoked by Vicsek and co-workers to mark a crossover 
from discontinuous to "continuous" transitions. 

In a recent letter (i^l, Aldana et al. study order- 
disorder phase transitions in random network models and 
show that the nature of these transitions may change 
with the way noise is implemented in the dynamics (they 
consider the angular and vectorial noises defined here). 
Arguing that these networks are limiting cases of Vicsek- 
like models, they claim that the conclusions reached for 
the networks carry over to the transition to collective mo- 
tion of the VM-like systems. They conclude in particular 
that in the case of "angular" noise the transition to collec- 
tive motion is continuous. We agree with the analysis of 
the network models, but the claim that they are relevant 
as limits of Vicsek-like models is just wrong: the data pre- 
sented there (Figure 1 of [40,] ) to substantiate this claim is 
contradicted by our Figures [5^,0 (see also [13]) obtained 
at larger system sizes. Again, for large-enough system 
sizes, the transition is indeed discontinuous. Thus, at 
best, the network models of Aldana et al. constitute a 
singular — > 00 limit of Vicsek-like models. 



(a) (b) 




FIG. 11: (color online) Typical snapshots in the ordered 
phase. Points represent the position of individual particles 
and the red arrow points along the global direction of mo- 
tion. (a)-(c): Angular noise, p = 1/2, vo — 0.5, i] = 0.3 and 
increasing system sizes - respectively L = 64, L = 256 and 
L — 1024. Sharp bands can only be observed if L is larger 
than the typical band width w. (d) Vectorial noise: p — 1/2, 
170 ~ 0.5, 77 — 0.55 and L — 64: bands appear at relatively 
small system sizes for this type of noise. For clarity reasons, 
only a representative sample of 10000 particles is shown in 
(b) and (c). Boundary conditions are periodic. 



IV. NATURE OF THE ORDERED PHASE 

We now turn our attention to the ordered, symmetry- 
broken phase. In previous analytical studies, it has often 
been assumed that the density in the ordered phase is 
spatially homogeneous, albeit with possibly large fluc- 
tuations (see, e.g. This is indeed what has been re- 
ported in early numerical studies, in particular by Vicsek 
et al. 0. In the following, we show that this is not true 
in large enough systems, where, for a wide range of noise 
amplitudes near the transition point, density fluctuations 
lead to the formation of localized, travelling, high-density 
and high-order structures. At low enough noise strength, 
though, a spatially-homogeneous ordered phase is found, 
albeit with unusually strong density fluctuations. 



A. Traveling in bands 

Numerical simulations of the ordered phase dynamics 
iv < rjt), performed at large enough noise amplitudes, 
are characterized by the emergence of high-density mov- 
ing bands {d = 2) or sheets (d = 3). Typical examples are 
given in Figs. llllI12l These moving structures appear for 
large-enough systems after some transient. They extend 
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FIG. 12: (color online) Same as Fig. [Tl] but in different ge- 
ometries/boundary conditions or space dimensions, (a)-(b) 
vectorial noise (77 = 0.325, p = 1/8, and vo = 0.5), boundary 
conditions are periodic along the y (vertical) axis and reflect- 
ing in X. (a): a long single band travels along the periodic 
direction, (b): the domain size along the periodic direction is 
too small to accomodate bands, and a single band bouncing 
back and forth along the non-periodic direction is observed, 
(c) Angular noise, repulsive force, and periodic boundary con- 
ditions (p = 2, ?7 = 0.23 and vo — 0.3). (d) High-density sheet 
traveling in a three-dimensional box with periodic boundary 
conditions (angular noise with amplitude 77 — 0.355, p = 1/2 
and Vo = 0.5). 



(a) (b) 




FIG. 13: (color online) (a) Typical density (black line) and 
order parameter (dashed red line) profiles for bands in two 
dimensions (vectorial noise, p = 2, — 0.6 and vo = 0.5). 
(b) Tail of the density profile shown in (d) (black line) 
and its fit (blue dashed line) by the formula: p±{x^^,t) ~ 
ao + ai{t) exp(— a:||/u;), with w « 6.3 (lin-log scales), (c-d) 
Traveling sheet in three dimensions (angular noise, p — 1/2, 
■q = 0.355, and vo = 0.5). (c) projection of particle positions 
on a plane containing the global direction of motion (marked 
by red arrow), (d) Density (black line) and order parameter 
(dashed red line) profiles along the direction of motion 



transversally with respect to the mean direction of mo- 
tion, and have a center of mass velocity close to vq. While 
particles inside bands are ordered and, in the asymptotic 
regime, move coherently with the global mean velocity, 
particles lying outside bands — in low density regions — 
are not ordered and perform random walks. 

As shown in Fig. [TTti-c for angular noise dynamics H]), 
there exists a typical system size Lb, below which the 
bands or sheets cannot be observed. Numerical simula- 
tions indicate that Lb depends only weakly on the noise 
amplitude and is of the same order of magnitude as the 
crossover size marking the appearance of the discontinu- 
ous character of the transition: Lb ~ L*. It is therefore 
numerically easier to observe bands in the ordered phase 
of vectorial noise dynamics ([3]), as in FigfTTH. 

Bands may be observed asymptotically without and 
with a repulsive interaction (Fig[T2t) and for both kinds 
of noise. They appear for various choices of boundary 
conditions (see for instance Fig. [T2k-b. where reflecting 
boundary conditions have been employed), which may 
play a role in determining the symmetry-broken mean 
direction of motion. For instance, bands travelling paral- 
lel to one of the axis are favoured when periodic boundary 
conditions are employed in a rectangular box (they repre- 
sent the simplest way in which an extended structure can 
wrap around a torus, and are thus reached more easily 



from disordered initial conditions), but bands travelling 
in other directions may also appear, albeit with a smaller 
probability. 

Bands can be described quantitatively through local 
quantities, such as the local density pi{x,t), measured 
inside a domain V{x) centered around x, and the local 
order parameter 

<fe{x,t) = — \{vl{t))p.^v(S)\ ■ (11) 
Vo 

Further averaging these local quantities perpendicularly 
to the mean velocity ([7]) one has the density profile 
p±{x^\,t) — {pi{x,t))± and the order parameter profile 
(p±{x\\,t) = {(pi{x,t)) ±, where x\\ indicates the longitu- 
dinal direction w.r.t. mean velocity. Bands are charac- 
terized by a sharp kink in both the density and the order 
parameter profiles (see Fig. [T5b and [T5b-d). They are 
typically asymmetric, as it can be expected for moving 
structures, with a rather sharp front edge, a well defined 
mid-height width w - which typically is of the same order 
as Lf, - and an exponentially decaying tail with a char- 
acteristic decay length of the order of the w (Fig. [T^) . 

Large systems may accomodate several bands at the 
same time, typically all moving in the same direction (see 
for instance Fig. [TTb and the density profile in Fig. [Hb). 
However, they do not form well-defined wave trains, but 
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rather a collection of solitary objects, as hinted by the 
following numerical experiments. 
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FIG. 14: (color online) Emergence of high-density high- 
order traveling bands (d — 2) from a spatially-homogeneous 
(uniformely distributed random positions) initial condition 
with all particle velocities oriented along the major axis of 
a 196 X 1960 domain with periodic boundary conditions. 
Vectorial noise of amplitude rj = 0.6, density p = 2, and 
Wo = 0.5. (a): space-time plot of the density profile (time 
is running from left to right from t = Q %o t = 12000, 
while the longitudinal direction is represented in ordinates. 
Colour scale from blue (low values) to red (high values).). 

(b) : same as (a) but at later times (from t = 148000 to 
t = 160000). (c): spatial Fourier power spectrum S of an 
early density profile (t = 12000). (d): early time evolution 
of selected Fourier modes (chosen between spectrum peaks of 

(c) ), k = 10, 23, 28, 33, 41, 76, 121 (inset: average over 50 dif- 
ferent runs), (e): density profile at a late time {t — 160000, 
final configuration of (b)). (f): same as (c) but for the late 
time density profile of (e) . The quality of this figure has been 
reduced to meet the ArXiv size constraint. 



We investigated the instability of the density- 
homogeneous, ordered state in a series of numerical sim- 
ulations starting from particles uniformly distributed in 
space but strictly oriented along the major axis in a large 
rectangular domain. Figure [T4k .b show space-time plots 
of the density profile: initially flat, it develops structures 
with no well-defined wavelength (Fig. [T4b). Density fluc- 
tuations destroy the initially-ordered state in a rather 
unusual way: a dynamical Fourier analysis of the density 
profile show a weakly-peaked, wide band of wavelengths 
growing subexponentially (Fig. I14d ). This is at odds 
with a finite-wavelength supercritical instability, which 
would lead to a wavetrain of traveling bands. Further- 
more, the asymptotic (late time) power spectra of the 
density profiles are not peaked around a single frequency 
either, but rather broadly distributed over a large range 
of wavenumbers (Fig. I14f ). In the asymptotic regime, 
bands are extremely long-lived metastable (or possibly 
stable) objects, which are never equally-spaced (a typi- 
cal late-time configuration is shown in Fig. [Hb). 

To summarize, the emerging band or sheet structure in 
the asymptotic regime is not a regular wave train char- 
acterized by a single wavelength, but rather a collection 
of irregularly-spaced localized traveling objects, probably 
weakly interacting through their exponentially decaying 
tails. 



B. Low-noise regime and giant density fluctuations 

As the noise amplitude is decreased away from the 
transition point, bands are less sharp, and eventually dis- 
appear, giving way to an ordered state characterized by 
an homogeneous local order parameter and large fluctu- 
ations of the local density. 

A quantitative measure of the presence, in the ordered 
phase, of structures spanning the dimension transverse 
to the mean motion (i.e. bands or sheets) is provided by 
the variances of the density and order parameter profiles: 

Aplit) = ((pi(a:||,i)-(pi(x[|,i))||)2)|| 

A^i(t) = ((^^(x||,t)-(^^(:r:||,t))||)2)|| ^'') 

where (•)|| indicates the average of the profile in the longi- 
tudinal direction with respect to mean velocity. Indeed, 
these profile widths vanish in the infinite-size limit except 
if band/sheet structures are present. 

In Fig [T5k.b. we plot these profile widths averaged 
over time as a function of noise amplitude. Both quan- 
tities present a maximum close to the transition point 
in the ordered phase, and drop drastically as soon as 
the disordered phase is entered. Lowering the noise 
away from the transition point, these profiles decrease 
steadily: bands/sheets stand less sharply out of the dis- 
ordered background (Fig [T5b). At some point (77 « 0.3 
for the parameters values considered in FigfTSb.b). bands 
rather abruptly disappear and are no more well-defined 
transversal objects. It is difficult to define this point ac- 
curately, but it is clear that for lower noise intensities the 
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FIG. 15: (color online) (a-b) Time- averaged profile width for 
both density (a) and order parameter (b) as a function of noise 
amplitude in the ordered phase (angular noise, 1024 x 256 
domain, global motion along the major axis, p = 1, and no = 
0.5). The dashed vertical blue line marks the order-disorder 
transition, (c-d) Typical instantaneous profiles along the long 
dimension of the system described in (a-b) for intermediate 
noise value ((b) r\ = 0.4) and in the bandless regime ((c) 
r] = 0.15). 
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FIG. 16; (color online) Giant density fluctuation and trans- 
verse superdiffusion in the bandless ordered phase. (a) 
anomalous density fluctuations (see text): An scales approx- 
imately like n" * (the dashed line has slope 0.8) both in two 
dimensions (black circles, L = 256, p = 2, uo = 0.5, angu- 
lar noise amplitude 77 = 0.25) and in three dimensions (red 
triangles, L — 64, p = 1/2 vo = 0.5, vectorial noise ampli- 
tude rj = 0.1, values shifted for clarity), (b) average doubling 
time T2 of the transverse (w.r.t. mean velocity) interparticle 
distance S±. Black circles: ordered bandless regime (p — 4, 
angular noise amplitude i] = 0.2 in a rectangular box of size 
1024 X 256) . The black dashed line marks the expected growth 



phase (p — 4, angular noise amplitude r] = 1, L — 512). 
dot dashed red line shows normal diffusive behavior: T2 ' 



The 



local order parameter is strongly homogeneous in space. 
Nevertheless, fluctuations in the density field are strong 
(Fig[T5tl). but can no more give rise to (meta)stable long- 
lived transverse structures. 

Density fluctuations in the bandless regime are in fact 
anomalously strong: measuring number fluctuations in 
sub-systems of linear size £, wc find that their root mean- 
square An does not scale like the square root of n = pi"^, 
the mean number of particles they contain; rather we 
find a{n) cx n" with a 0.8 both in d = 2 and in 
d = 3 (Fig. [TBk). This is reminiscent of the recent dis- 
covery of "giant density fluctuations" in active nemat- 
ics [3, HH, [54I ■ However, the theoretical argument which 
initially predicted such fluctuations [s^ cannot be in- 
voked directly in the present case. (Indeed, the above 
value of a, although needing to be refined, does not seem 
to be compatible with the prediction Oi — ^ + ^ niade 
in [5^.) More work is needed to fully understand un- 
der which circumstances the coupling between density 
and order in systems of "active" , self-propelled particles 
gives rise to such anomalous density fiuctuations. 

C. Transverse superdiffusion 

According to Toner and Tu predictions [1, [H, [13] , the 
dynamics of the symmetry-broken ordered phase of polar 
active particles should be characterized by a superdiffu- 
sive mean square displacement 

Ax^ = V((a;i(i)-a;±(0))2), (13) 



in the direction(s) transversal to mean velocity. In par- 
ticular, in d = 2 one has [s^ 

Axl - r (14) 

with ly = 4/3. While this analytical result has been suc- 
cesfully tested by numerical simulations of models with 
cohesive interactions [sH, [38|, numerical simulations in 
models without cohesion present substantial difficulties, 
mainly due to the presence of continuously merging and 
splitting sub-clusters of particles moving coherently (as 
discussed in Section [IV D[) . As a consequence, an ensem- 
ble of test particles in a cohesionless model is exposed 
to different "transport" regimes (w.r.t. center of mass 
motion) which are not well separated in time. When the 
mean displacement is averaged at fixed time, this tends 
to mask the transverse superdiffusion. 

To overcome this problem, we chose to follow [s^] and 
to measure T2, the average time taken by two particles 
to double their transverse separation distance S±. From 
Eq. one immediately has 

r2 - (15) 

with 2lv — 3/2 va. d — 2. In order to easily separate the 
transverse from the parallel component, we considered an 
ordered system in a large rectangular domain with peri- 
odic boundary conditions and the mean velocity initially 
oriented along the long side. The mean direction of mo- 
tion then stays oriented along this major axis, so that we 
can identify the transverse direction with the minor axis. 
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Furthermore, a high density and a small (angular) noise 
amplitude (corresponding to the bandless regime) have 
been chosen to avoid the appearance of large, locally dis- 
ordered patches. 

Our results (Fig. [T6b) confirm the prediction of Toner 
and Tu: transverse superdiffusion holds at low-enough 
noise, while normal diffusion is observed in the disor- 
dered, high-noise phase. Note that the systematic devia- 
tion appearing in our data at some large scale is induced 
by large fluctuations in the orientation of the global mean 
velocity during our numerical simulations (not shown 
here). 

We take the opportunity of this discussion to come 
back to the superdiffusive behavior of particles observed 
in the transition region[2lj. There, sub-clusters emerge 
propagate ballistically and isotropically due to the ab- 
sence of a well-established global order. Particles trajec- 
tories consist in "ballistic flights" , occurring when a par- 
ticle is caught in one of these coherently moving clusters, 
alternated with ordinary diffusion in disordered regions. 
The mean square displacement of particles exhibits the 
scaling Ax^ = {\x{t) - ^(0))^)^ cx t^/s [2l[. In view of 
our current understanding of the discontinuous nature of 
the transition, we now tend to believe that this isotropic 
superdiffusion is probably not asymptotic. 



D. Internal structure of the ordered region 

We now turn our attention to the internal structure of 
the ordered regimes. As we noted in the previous sec- 
tion, these regimes do not consist of a single cluster of 
interacting particles moving coherently. Even in the case 
where high density bands/sheets are present, these are 
in fact dynamical objects made of splitting and merg- 
ing clusters. Note that for the models considered here, 
clusters are unambiguously defined thanks to the strictly 
finite interaction range tq. 

As noticed first by Aldana and Huepe [7,], clusters of 
size n are distributed algebraically in the ordered region, 
i.e. P{n) ^ n~^. But a closer look reveals that the ex- 
ponent /X characterizing the distribution of cluster sizes 
changes with the distance to the transition point. For 
noise intensities not too far from the threshold, when 
bands are observed, we find /i values larger than 2, 
whereas /i < 2 in the bandless regimes present at low 
noise intensities (Fig. [TTk.b). 

Thus, bands are truly complex, non-trivial structures 
emerging out of the transverse dynamics of clusters with 
a well defined mean size (since /i > 2). It is only in 
the bandless regime that one can speak, as Aldana and 
Huepe, of "strong intermittency" . We note in passing 
that the parameter values they considered correspond in 
fact to a case where bands are easily observed (at larger 
sizes than those considered in [7|). Thus, clusters do 
have a well-defined mean size in their case. Consequently, 
the probability distribution P{ip) of the order parameter 
ip does not show the behavior reported in Fig. 1 of [7| 
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FIG. 17: (a-b) Cluster size distributions (arbitrary units) for 
domain sizes L = 32 (black), 128 (red) , and 512 (green) from 
left to right {d = 2, p = 2, vo = 0.5, angular noise), (c-d) 
Probability distribution functions of the order parameter ip 
(arbitrary units) for the same parameters and system sizes as 
in (a-b). (The most peaked distributions are for the largest 
size L — 512.) Left panels (a,c): r] = 0.1, bandless regime; 
Right panels (b,d), regime with bands at = 0.4. 



as soon as the system size is large enough. Whether 
in the band/sheet regime or not, P{ip) shows essentially 
Gaussian tails, is strongly peaked around its mean, and 
its variance decreases with system size fFig. [TTb.d). 

Although the picture of intermittent bursts between 
"laminar" intervals proposed by Aldana and Huepe has 
thus to be abandoned, the anomalous density fluctua- 
tions reported in the previous section are probably tanta- 
mount to the strong intermittency of cluster dynamics in 
the bandless regime. Again, these phenomeiia, reported 
also in the context of "active" nematics [3, [5l|, [s^l , de- 
serve further investigation. 



E. Phase ordering 

The ordered regimes presented above are the result of 
some transient evolution. In particular, the bands/sheets 
are the typical asymptotic structures appearing in finite 
domains with appropriate boundary conditions. In an in- 
finite system, the phase ordering process is, on the other 
hand, infinite, and worth studying for its own sake. 

Numerically, we have chosen to start from highly disor- 
dered initial conditions which have a homogeneous den- 
sity and vanishing local order parameter. In practice, 
we quench a system "thermalized" at strong noise to a 
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FIG. 18: (color online) Phase-ordering from disordered initial conditions (d = 2, angular noise amplitude ry = 0.08, p = 1/8, 
vo = 0.5, system size L = 4096). Snapshots of the density field coarse-grained on a scale ^ = 8 at times t — 160, 320, and 640 
from left to right. 
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FIG. 19: (color online) Phase ordering as in Fig.[T8](I/ = 4096 
p — 1/2, vq — 0.5). (a) two-point density correlation function 
C{r,t) — {pe{x,t)pi{x + r,t))s (coarse grained over a scale 
£ = 4) as a function of distance r = jr| at different time steps: 
from left to right t = 50, 100, 200, 400, 800, 1600. Noise ampli- 
tude is r; = 0.25, data have been further averaged over « 40 
different realizations. Inset: log scales reveal the interme- 
diate near-algebraic decay and the quasi-exponential cut-off. 
(b) Lengthscale ^, estimated from the exponential cut-off po- 
sitions, as a function of time. Empty black circles: 77 — 0.25 
as in (a) (i.e. regime in which bands are observed asymp- 
totically). Red full triangles: = 0.1 (i.e. in the bandless 
regime). The dashed black line marks linear growth. 



smaller, subcritical, 77 value. Typical snapshots show the 
emergence of structures whose typical scale seems to in- 
crease fast (Fig. \TE\\ . During this domain growth, we 
monitor the two-point spatial correlation function of both 
the density and velocity fields. These fields are defined 
by a coarse-graining over a small lengthscale £ (typically 
4). These correlation functions have an unusual shape 
(Fig. [TOb): after some rather fast initial decay, they dis- 
play an algebraic behavior whose effective exponent de- 
creases with time, and finally display a near-exponential 
cut-off. As a result, they cannot be easily collapsed on 
a single curve using a simple, unique, rescaling length- 



scale. Nevertheless, using the late exponential cut-ofF, a 
correlation length ^ can be extracted. Such a lengthscale 
^ grows roughly linearly with time (Fig. (T^b). Quali- 
tatively similar results are obtained whether the noise 
strength is in the range where bands/sheets appear in 
finite boxes or not. 

We note that the above growth law is reminiscent of 
that of the so-called model H of the classification of 
Halpcrin and Hohcnberg [s^ . Since this model describes, 
in principle, the phase separation in a viscous binary 
fluid, the fast growth observed could thus be linked to 
the hydrodynamic modes expected in any continuous de- 
scription of Vicsek-like models [1, HH . 



V. GENERAL DISCUSSION AND OUTLOOK 
A. Summary of main results 

We now summarize our main results before discussing 
them at a somewhat more general level. 

We have provided ample evidence that the onset of 
collective motion in Vicsek-style models is a discontinu- 
ous (first-order) phase transition, with all expected hall- 
marks, in agreement with We have made the (nu- 
merical) effort of showing this in the limits of small and 
large velocity and/or density. 

We have shown that the ordered phase is divided in 
two regions: near the transition and down to rather 
low noise intensities, solitary structures spanning the di- 
rections transverse to the global collective velocity (the 
bands or sheets) appear, leading to an inhomogeneous 
density field. For weaker noise, on the other hand, no 
such structures appear, but strong, anomalous density 
fluctuations exist and particles undergoes superdiffusive 
motion transverse to the mean velocity direction. 

Finally, we have reported a linear growth (with time) 
of ordered domains when a disordered conflguration is 
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quenched in the ordered phase. This fast growth can 
probably be Unked to the expected emergence of long 
wavelength hydrodynamic modes in the ordered phase of 
active polar particles models. 



B. On the role of bands/sheets 



The high-density high-order traveling bands or sheets 
described here appear central to our main findings. They 
seem to be intimately linked to the discontinuous charac- 
ter of the transition which can, to some extent, be consid- 
ered as the stability limit of these objects. In the range 
of noise values where they are observed, the anomalous 
density fluctuations present at lower noise intensities are 
suppressed. 

One may then wonder about the universality of these 
objects. Simple variants of the Vicsek-style models stud- 
ied here (e.g. restricting the interactions to binary ones 
involving only the nearest neighbor) do exhibit bands 
and sheets [ll|. Moreover, the continuous deterministic 
description derived by Bertin et al. 'QJ^ does possess lo- 
calized, propagating solitary solutions rather similar to 
bands [53|. Although the stability of these solutions 
need to be further investigated, these results indicate 
that these objects are robust and that their existence is 
guaranteed beyond microscopic "details". However, the 
emergence of regular, stable, bands and sheets is obvi- 
ously conditioned to the shape and the boundary con- 
ditions of the domain in which the particles are allowed 
to move. In rectangular domains with at least one pe- 
riodic direction, these objects can form, span across the 
whole domain, and move. But in, say, a circular domain 
with reflecting boundary conditions, they cannot develop 
freely, being repeatedly "frustrated". Nevertheless, sim- 
ulations performed in such a geometry indicate that the 
transition remains discontinuous, with the ordered phase 
consisting of one or several dense packets traveling along 
the circular boundary. Note, though, that these packets 
intermittently emit elongated structures (bands) travel- 
ing towards the interior of the disk before colliding on the 
boundary. To sum up, bands appear as the "natural" ob- 
jects in the transition region, but they may be prevented 
by the boundaries to develop into full-size straight ob- 
jects. 

At any rate, time-series of the order parameter such 
as the one presented in Fig. [5] clearly show that the 
transition is discontinuous irrespective of the geometry 
and boundaries of the domain, and thus of whether 
bands/sheets can develop into stable regular structures or 
not: the sudden, abrupt, jumps from the disordered state 
to some ordered structure are tantamount to a nucleation 
phenomenon characteristic of a discontinuous transition. 
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FIG. 20: (color online) Scatter plots of local order parameter 
vs local density in the ordered phase (angular noise, p — 2, 
vo — 0.5 in a domain of size 1024 x 256, (same parameters 
as in Fig. 115b). The local quantities were measured in boxes 
of linear size ^ = 8. Here, the local order is represented by 
the angle between the orientation Qiocai of the local order pa- 
rameter and the global direction of motion Q global- The black 
solid lines are running averages of the scatter plots. The red 
solid lines indicate the global density p = 2 (and thus marks 
the percolation threshold in a two-dimensional square lattice). 
(a,b): Tj — 0.1 and 0.2: in the bandless regime, the ordered 
plateau starts below p = 2, i.e. ordered regions percolate, (c) 
approximately at the limit of existence of bands: the start of 
the plateau is near p = 2. (d) at higher noise amplitude in 
the presence of bands. 



C. A speculative picture 

We would now like to offer the following speculative 
general picture: the key feature of the Vicsek-like models 
studied here — as well as of other models for active me- 
dia made of self-propelled particles [ll|, [H, [H, [H, — 
is the coupling between density and order. Particles are 
forced to move, and, since they carry informations about 
the order, advection, density fluctuations and order are 
intimately linked. High-density means strong local order 
(if noise is low enough) because the many particles in a 
given neighborhood will adopt roughly the same orienta- 
tion. The reverse is also true: in a highly-ordered region, 
particles will remain together for a long time and thus 
will sweep many other particles leading to a denser and 
denser group. 

At a given noise level, one can thus relate, in the spirit 
of some local equilibrium hypothesis, local density to lo- 
cal order. In practice, such an "equation of state" ap- 
proach can be justified by looking, e.g., at a scatter plot 
of local order parameter vs local density. Fig. [201 re- 
veals that, in the ordered bandless regime, such a scat- 
ter plot is characterized by a plateau over a large range 
of local density values corresponding to order, followed, 
below some crossover density, by more disordered local 
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patches. The regions in space where local density is be- 
low this crossover level do not percolate in the band- 
less regime, and order can be maintained very steadily 
in the whole domain (this is corroborated by the fact 
that in spite of the large, anomalous density fluctua- 
tions, the order parameter field is, on the other hand, 
rather constant, see Fig. [T5b). The noise intensity at 
which bands emerge roughly corresponds to the value 
where the low-density disordered regions percolate. The 
reamining disconnected, dense patches then eventually 
self-organize into bands/sheets. The emergence of these 
elongated structures is rather natural: moving packets 
elongate spontaneously because they collect many par- 
ticles; superdiffusion in the directions perpendicular to 
the mean motion endow these nascent bands/sheets with 
some "rigidity". At still stronger noise, the bands/sheets 
are destroyed and global order disappears. 

The above features are at the root of the approach 
by Toner and Tu ^8,, ,22n ^] ■ Their predictions of strong 
density fluctuations, transverse superdiffusion, and pecu- 
liar sound propagation properties are correct as long as 
bands/sheets do not exist, i.e. for not too strong noise in- 
tensities. This is indeed in agreement with their assump- 
tion that the density field is statistically homogeneous 
in the frame moving at the global velocity (albeit with 
strong fluctuations), which is only true in the bandless 
regime. 



D. Outlook 

The results presented here are almost entirely numer- 
ical. Although they were obtained with care, they need 
to be ultimately backed up by more analytical results. 
A first step is the derivation of a continuous description 
in terms of a density and a velocity field (or some com- 
bination of the two), which would allow to go beyond 
microscopic "details". In that respect the deterministic 
equation derived by Bertin et al. from a Boltzmann de- 
scription in the dilute limit [ll[ is encouraging. However, 
one may suspect that intrinsic fluctuations are crucial in 
the systems considered here if only because some of the 
effective noise terms will be multiplicative in the density. 
A mesoscopic, stochastic equation description is thus a 
priori preferable. This is especially true in view of the 
"giant" anomalous fluctuations present in the bandless 
ordered phase. These fluctuations clearly deserve further 
investigation, all the more so as they seem to be generic 
features of active particle models [13, [IH . 

Ongoing work is devoted to both these general issues. 
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FIG. 21; (color online) Test of the fluctuation-dissipation re- 
lation on the vectorial model with repulsive force (p — 2, 
vq = 0.3, L = 128). (a): susceptibility vs time at reduced 
noise amplitude e = 1 — r]/rjt = 0.005, \h\ — 10~^ (dashed red 
line) and \h\ — 10^'^ (plain black line), (b); susceptibility vs. 
correlation at \h\ = 10"^ and e = 0.08, 0.26 and 0.44 from 
top to bottom, (c): effective temperature vs noise amplitude. 
The vertical dashed line marks the transition point. 



edged. We thank A. Vulpiani and M. Cencini for intro- 
ducing us to the method outlined in Ref. [sS ]. 



APPENDIX A: FLUCTUATION-DISSIPATION 
RELATION 

In 0, Vicsek et al. also studied the validity of the 
fluctuation-dissipation theorem and concluded, from nu- 
merical analysis, that it is violated. Here we approach 
this question again, using the "dynamical" approach put 
forward in [60], rather than the equilibrium used in 
The fluctuation-dissipation relation is expressed as: 



R{t - to) 



1 dC{t-to) 



off 



dt 



(Al) 



where R is the response, C the associated correlation 
function, and T^s is some "effective temperature" . 

In both cases, at any rate, the key point is to investi- 
gate the effect of an external field on the ordering process. 
Since here one cannot rely on any Hamiltonian structure, 
the external field remains somewhat arbitrary, as it can- 
not be unambiguously defined as the conjugate variable 
of the order parameter. 

In the field h was directly confronted to the local 
average velocity, and the governing equation was replaced 
by: 



V,{t + M) = Vo {TZrjOd) 



(A2) 
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This manner of introducing h leads to an effective in- 
tensity which depends on the local ordering: h is com- 
paratively stronger in disordered regions (an important 
effect at the early stages of ordering) than in ordered 
regions. This could in fact prevent the necessary linear 
regime to occur, even at very low field values (see, for 
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instance, Fig. 6 of Q). A mean- field analysis has con- 
firmed this view, showing a logarithmic variation of the 
response with the field intensity [61] . 

To bypass this problem, we have preferred to use the 
following equation: 



vl{t + At) ^voi9 



jeVi 



h + 



(A3) 



where the vectorial noise was chosen because, as shown 
above, it leads more easily to the asymptotic regime. The 
effective intensity of the field is now proportional to the 
local order. 

Two (scalar) response functions can be defined in our 
problem: the longitudinal response along the field di- 
rection, and the transverse response R±. We consider the 
former. In practice, we quenched, at time Iq, a strong- 
noise, highly-disordered system ip(to) ~ to a smaller 
noise value and started applying the constant, homo- 
geneous field h immediately. We then followed the re- 
sponse of the system by monitoring the growth of the 
order parameter. We measured the susceptibility X|| 
which is nothing but the integrated response function: 
X||(t,to) = /(* R\\{t,t')dt' . In practice, we have: 



x\\{t-h) = -^0{t)-m 
\h\ 



(A4) 



In a well-behaved system, the susceptibility should be in- 
dependent of the amplitude of the field, at least at small 
enough values ("linear" regime). This is what we ob- 
served, as shown in Fig. [^Tb . 

Correspondingly, the correlation function is defined as: 



C{t ~ to) 



1 



[vi{to) ■ ni.{t))t . 



(A5) 



The fluctuation-dissipation relation 
written in its integrated form: 



Jeff 



Til can then be 



(A6) 



C(t-to)) 



In Fig. [2Tb. we show that x\\ and C are related linearly in 
time, conflrming the validity of this relation and allowing 
an estimation of Tcs- 

This well-defined — although not uniquely defined — 
effective temperature varies as expected in parameter 
space. In particular, it increases systematically with the 
noise strength rj (Fig. \21k). although this variation is not 
linear. Note also, that, intriguingly, there is a small jump 
of Toff at the noise value corresponding to the transition 
in this case {rjt ~ 0.55). 
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